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Chaos, Determinacy and Fractals in Active-Sterile Neutrino Oscillations 

in the Early Universe 

Kevork N. Abazajiai^j and Prateek Agrawajj 
Department of Physics, University of Maryland, College Park, MD 20742-4111 

The possibility of light sterile neutrinos allows for the resonant production of lepton number in 
the early universe through matter-affected neutrino mixing. For a given a mixing of the active 
and sterile neutrino states it has been found that the lepton number generation process is chaotic 
and strongly oscillatory. We undertake a new study of this process' sensitivity to initial conditions 
through the quantum rate equations. We confirm the chaoticity of the process in this solution, and 
moreover find that the resultant lepton number and the sign of the asymmetry produces a fractal 
in the parameter space of mass, mixing angle and initial baryon number. This has implications for 
future searches for sterile neutrinos, where arbitrary high sensitivity could not be determinate in 
forecasting the lepton number of the universe. 

PACS numbers: 14.60.Pq; 14.60.St; 26.35. +c; 95.30.-k 



I. INTRODUCTION 



Neutrino physics is entering an era of increasing pre- 
cision as initial evidence for neutrino oscillations is lead- 
ing to confirmation and refinement of the neutrino mass 
and mixing parameters. The solar neutrino deficit 
is consistent with a matter-affected Mikheyev-Smirnov- 
Wolfenstein (MSW) transformation of electron neutrinos, 
z/e pLj, and was confirmed in the CPT conjugate channel 
by reactor electron antineutrino disappearance at Kam- 
LAND [5]. Atmospheric neutrino disappearance due to 
apparent neutrino oscillations between what is predomi- 
nantly i/f^ and Vr and their CPT conjugates has been con- 
firmed and is being determined precisely by long-baseline 
neutrino experiments '3', '5 . 

There were also indications of neutrino oscillations gov- 
erned by a large mass scale from the short-baseline Liquid 
Scintillator Neutrino Detector (LSND) "W that would re- 
quire a third mass scale difference and therefore a fourth 
neutrino or more [5] . This neutrino could not participate 
in weak interactions due to the Z^ width [7], and there- 
fore is dubbed "sterile." Sterile neutrinos are prevalent in 
most models of neutrino mass generation mechanisms [S] , 
and could be light enough to participate in neutrino oscil- 
lations. The four neutrino oscillation hypothesis for the 
LSND result is in conflict with new results from the from 
the short-baseline MiniBOONE experiment |9] . The pa- 
rameter space of interest for the LSND and MiniBOONE 
results are a small fraction of that of interest for lepton 
number generation, and do not exclude the possibility of 
this mechanism. Moreover, such light sterile neutrinos in 
the parameter space of interest for lepton number gener- 
ation, and outside of the MiniBooNE bounds are favored 
to produce the electron fraction required by successful 
r-process nucleosynthesis in Type-II supernovae flOl lllj . 
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See-saw type mass generation mechanisms could readily 
include such light sterile neutrinos [12]. 

The effects of sterile neutrinos in the early universe 
have been under study for some time. Early work stud- 
ied the effect of the partial or complete thermalization 
of the sterile neutrino through oscillations, which would 
alter the energy density and expansion rate of the uni- 
verse prior to primordial nucleosynthesis, and therefore 
the production of the light elements [131 US HSl US] • The 
primordial helium yield is very sensitive to the expan- 
sion rate, and has been used to place bounds on four 
neutrino models with appreciable mixing between active 
and sterile neutrinos of the kind required for the LSND 
indications [TT] [18] . It was also discovered that an MSW 
resonance can cause the generation of a large neutrino 
asymmetry that could modify primordial nucleosynthe- 
sis even with very small mixing angles between the active 
and sterile neutrinos [19] . This was confirmed in Ref. [20] , 
and further found rapid oscillations in the asymmetry 
generation, extreme sensitivity to initial conditions, and 
chaoticity in certain regions of parameter space of mass 
and mixing angle. 

A key question regarding the chaoticity of the lepton 
generation mechanism is the determinacy of the sign of 
the asymmetry given a mixing angle and mass-splitting 
for the active-sterile neutrino pair. That is, if a mass 
and mixing angle were experimentally measured for the 
pair, could the sign and magnitude of the cosmological 
neutrino asymmetry be predicted? One way of assess- 
ing the predictability of the final lepton number given a 
potentially experimentally determined mass and mixing 
angle is the level of information loss in the lepton number 
generation process. The level of information loss in the 
chaotic process was explored in Ref. [21], which found 
appreciable loss of initial conditions due to chaoticity in 
small regions of parameter space. The sensitivity of the 
sign and amplitude of the asymmetry to the mixing pa- 
rameter space was found to be "complex" in Ref. [22] . 

All of the work discussed thus far was done through 
the quantum rate equations (QRE's), where the neutrino 



system is modeled to follow the average neutrino momen- 
tum (p) « 3.1 ST. The strong momentum dependence of 
the resonant process of lepton number generation sug- 
gests that an accurate solution should take into account 
the full spectrum in the solution of the equations through 
the so-called Quantum Kinetic Equations (QKE's). The 
narrowness of the resonance requires exceptionally fine 
momentum binning and therefore following the exponen- 
tial growth in a very large number of momenta modes. 
Several authors found limitations to solving the equa- 
tions numerically [531 [H]. For example, Ref. [23] found 
oscillatory behavior in a range of parameter space with 
the QKE's, but could not conclude that it was physi- 
cal and not numerical. The QKE's were solved through 
the oscillatory phase by Ref. "55], who used a new bin- 
ning procedure following the resonance to simulaneously 
have fine binning yet a manageable number of momenta 
modes. They found that the chaotic nature of the gener- 
ation process remained. 

In this work, we explore the nature of the lepton num- 
ber generation process and its chaoticity. We confirm 
that the resonant generation is highly oscillatory in re- 
gions of the parameter space and extremely sensitive to 
the initial conditions. Moreover, we explore the depen- 
dence of the final sign of the asymmetry on the oscillation 
parameter space and find that the extreme sensitivity to 
the parameter space gives rise to a fractal structure in 
the final sign dependence on mass-splitting and mixing 
angle. This has significant implications for the determi- 
nacy of the cosmological lepton number in the case of 
discovering a small sterile neutrino mixing with the ac- 
tive neutrino sector. If the parameters are found to lie in 
a region of fractal nature, even arbitrarily high precision 
in the determination of the parameters will not allow for 
a prediction of the cosmological lepton number. 

We will briefly introduce our system in Sec. |ll] includ- 
ing the effects of our constant momentum approxima- 
tion. We go on in Sec. |III| to present the solution under 
the above approximation, where we show that the final 
sign of the asymmetry is not just chaotic, but exhibits a 
self-similarity like that of a fractal at very small scales. 
We compare our results with Ref. [H] , which shows that 
the resolution and numerical error are responsible for the 
seeming "white noise" structures in the parameter space. 
We conclude in Sec. IIVI 



An a-type neutrino asymmetry is defined as 



(2) 



where the photon number density n^ = 2({3)T^/tt'^. The 
effective total lepton number L^"-', useful when defining 
the evolution equations below, is given by 
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where rj = n\,ln^ is the baryon to photon ratio at the 
epoch of oscillations. 

The density matrix parametrization for an active- 
sterile oscillation can be chosen as 



Pr;3'(p) = ^[Po(p)/ + P(p).a], 



(4) 



where P(p) = Pj;(p)x 4- -P^y -I- P^^ is the polarization 
vector; g = crjjX + cryy-i-cr^z are the usual Pauli matrices. 
The evolution of P(p), Po(p)i P(p) and Po(p) is given by 
the following equations [25, , 



dV 

dPa 



= V{p) X P{p) - D{p)[P,{p)ic + Py{p)]y] 



dPo 
dt 



dt -nv)[^-\ip.iP)+p.m. 



(5) 



The equations for the anti-particles are given by the 
obvious substitutions, P -^ P, Pq -^ Pq, V(p) -^ V(p), 
feq{p) "> feq{p)- V(p) is obtained by replacing L^^^by 
-~L^'^\ The total collision rates of the fiavor neutrino 
and anti-neutrino are approximately equal, T(p) ~ r(p). 
These 8 equations form the quantum kinetic equations 
for active-sterile oscillation. 

The damping coefficient is given by 



D{p) = T{p)/2, 



(6) 



where T{p) is the total collision rate of the neutrino, given 
by 0.92G%T''p for the r-ncutrino. Further, D ~ D to a. 
good approximation J22J . 

The rotation vector V(p) has the following components 



II. THE QUANTUM KINETIC AND RATE 
EQUATIONS 

For specificity, we concentrate on oscillations between 
the r-neutrino (v-r) and a sterile neurino (vs)- The vac- 
uum mixing angle is defined by the following equation. 



cos (6*0) sin (^o) 
-sin (6*0) cos (6*0) 






(1) 



where v^ and v^, are mass eigenstates. Conventionally, 
cos{0) is chosen to be positive, and 5rn? = ml — m^. 
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Here 



is the photon equiHbrium number density, 



rij/^ and n^- are normalized to unity and feqip) is the 
Fermi-Dirac distribution with a chemical potential //t-. 

The chaotic nature of the solution is seen to not depend 
on 77, which was set to ^ 1 x 10"^*^ (though the sign of 
the asymmetry can be sensitive to this initial value as 
well). The well-known Fermi-Dirac distributions are 



foip) = 



1 



1 + exp(j5/r) 
^"''^P^ = l + eMp/T-HrV 



(8) 



where jir = ^^v^/T. To get the distribution for anti- 
neutrinos, we substitute jlr ^ i^t- The chemical poten- 
tials can be found using the neutrino asymmetry. For a 
system in thermal equilibrium, the asymmetry is given 

by 
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X dx 
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X dx 



(9) 



We can expand the above integral in (jIt — jlf) 

+ 6 (A? - Al) ln(2) + (A? - A?) • (10) 



We need another relation between JIt and jlf to solve 
for them. For temperatures greater than the chemical 
decoupling temperatures, pair creation and annihilation 
processes are fast enough to ensure chemical equihbrium. 
However, below these temperatures, neutrinos no longer 
remain in thermal equilibrium, though weak interactions 
are fast enough to thermalize them. As a consequence, 
the value of /i-r remains fixed at the decoupling temper- 
ature, but jlf continues to decrease. In summary 



/i^ + jlf = Q T> rj"^^ 
jlr = jlr{T = TJgJ T < Tj;^, 



(11) 



where TJ^^ ~ 3.5 MeV. These two equations above de- 
termine the chemical potentials to be used in the Fermi- 
Dirac distributions. 

We employ the approximation that all the neutrinos 
have the same momentum, viz. p = {p) ~ 3.15T for our 
analysis. The approximation inflates the rate at which 
neutrino asymmetry is created and therefore may be un- 
suitable for the thermal regions where the asymmetry 
growth is shown to be exponential |23j . Apart from this 
shortcoming, this approximation has the signiflcant ad- 
vantage of being numerically tenable, and is useful to ex- 
plore the general features of the system. We also carried 
out analysis using the full momentum dependent quan- 
tum kinetic equations (QKEs). The full system is too 
numerically intensive to run over the entire parameter 
range at high resolution to study the fractal structure. 
Moreover, our work with the full QKEs never reached 
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FIG. 1: Shown is the evolution of the lepton number over 
temperature for a particular case of Vr ^^ Vs mixing, 5m? = 
-0.01 eV, sin^ 29 = 2 x 10"^ 



the numerical stability which we required, and describe 
below. 

Further, following the work of Ref. |22j . we approxi- 
mate Pq to be a constant, set to 1. The asymmetry is 
given simply by 



L,. 



{Pz - Pz 



(12) 



In order to achieve numerical efficiency and stability, 
we rewrite {Pi, Pi} in the above equation as {P^,P^}, 
where P^ = Pi ± Pi. This separates out the large quan- 
tities and the small quantities in the differential equation 
P2] . Thus, the equations take the form 



(13) 



P^ = 


= -dp:^- 


- VoP+ - VlP- 


p^ - 


-- -Dp^- 


- V^P+ + VoP+ + VlP, 


Pt ~- 


= y^K 




p^ = 


= -DP- 


-VoP--V,P+ 



p:. 



-DP- - V^P, - VoP- + VlP+ 



P7 = V.P-. 

These are coupled non-linear differential equations, 
with no analytical solution. The system is seen to un- 
dergo resonance at Vo = 0, if Jm^ < 0. The asym- 
metry L^^ undergoes rapid oscillation at the resonance, 
finally settling down to a power-law-like growth. See 
Fig. [1] for an example of the oscillatory evolution. The 
magnitude of the final asymmetry is given by the so- 
lution oi Vl + Vo = 0, and is given by \L„^\ ~ 
11 {\5m^\/eY^) (MeV/T)"*. After resonance, L^^ = be- 
comes an unstable extremum, and the system is trapped 
into one of the two minima given above. A detailed dis- 
cussion of how the asymmetry evolves can be found in 
Ref. j22]. The final sign of the asymmetry, however, is 
not fixed. It is seen to be very sensitive to initial con- 
ditions, mixing parameters, even choice of step sizes and 
integrators [25]. 
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FIG. 2: Shown is the final sign of the lepton asymmetry as 
a function of mass splitting and mixing angle in the case of 
i/t ^^ i^a mixing for a wide range of parameter space (top), 
and a zoomed region (bottom). White indicates a positive 
sign and black negative. 



The differential equations (Eq. [5]) are expressed in 
terms of time-derivatives but from a cosmological point 
of view, we are interested to look at the solutions in terms 
of temperature. The time-temperature relation is given 
by dt/dT ~ -Mpi/5.UT^, Mpi being the Planck mass. 



III. RESULTS 

We plot the asymmetry as a function of temperature 
for a couple of parameter sets, and our results match up 
well with previous work |22j . When we scan the entire pa- 
rameter range shown in Fig. [2J we note some interesting 
facts. First, in our calculations, the band-structure seems 
to continue up to many more bands as compared to seen 
in Ref. |22j . Further, we can make out structures even 
beyond the primary band structure, which has previously 
been described as pure white noise. Our numerical rou- 
tines show that the number of oscillations increases as the 
value of ^0 increases, and hence the number of steps for 
the evolution of the differential equations increases. This 
requires that the absolute error per step in the differen- 
tial equation be very small in order for the final sign to 
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FIG. 3: Shown is the final sign of the lepton asymmetry as 
a function of mass splitting and mixing angle in a zoomed 
parameter range to show detail of the sign dependence. 



be numerically stable. Also, beyond the primary bands, 
the structure becomes very fine, and an analysis with a 
coarse resolution will just sample the fine structure as an 
apparent random sign distribution. 

To resolve the fine structure, we zoom into a part of the 
parameter space. We see that the region which shows no 
interesting features in entire range, shows definite struc- 
ture when we improve the resolution in Fig. |3] with fur- 
ther zoom in to the edge of the bands seen in figure. The 
edge is seen to have structure at very high zoom mag- 
nitudes. This behavior is replicated when we zoom into 
other parameter ranges. This suggests that the sign of 
asymmetry exhibits the characteristics of a fractal with 
respect to mixing parameters. 

We employed the Gnu Scientific Library extensively for 
our numerical analysis"'^. The system of differential equa- 
tions are extremely sensitive numerically, and therefore 
demand a very careful treatment. We use the adaptive 
Runge-Kutta-Fehlberg routine for our system. For the 
momentum averaged solutions, we use a very stringent 
error bound per step of evolution, so that the accumu- 
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lated error does not affect the final sign of asymmetry. 
We check that our results converge for a more stringent 
error bound. As ^o and —8m? increase, so does the num- 
ber of resonant oscillations. This in turn implies that 
the step size would decrease in order to follow the reso- 
nance correctly, giving us a larger number of steps overall. 
We kept our error bound near 1 x 10~^^, and we found 
that the error begins to dominate the solution around 
1 X 10^^''. This may be why we see significantly more 
structure in Fig. |2] than what was previously seen, in ad- 
dition to the high resolution sampling of the parameter 
space. 

To characterize the nature of the fractal solutions, we 
estimate its fractal dimension with the correlation inte- 
gral 



C(e) 



1 



N 



\\X,-Xj\\) 



(14) 



where xi are the positions of the positive sign lepton num- 
ber in parameter space and Q{x) is the Heaviside func- 
tion. For finite N the integral is the correlation sum. For 
fractal structures, at high N the sum tends to 



C(6) 



(15) 



where a is the correlation dimension, an approximation 
of the fractal dimension. For the range —0.003 cV^ < 
bm^ < -0.001 eV^, 10"^ < sm^{29) < 2 x 10"^ we find a 
power-law dependence of the correlation sum, and arrive 
at a convergence of the correlation dimension of approx- 
imately 1.05 for N > 10^. For larger mixing angles, we 
obtain higher correlation dimensions of 1.35. 



IV. CONCLUSION 

We have explored the generation of cosmological lep- 
ton number in the case of the presence of a light ster- 
ile neutrino state with small mixing with an active neu- 
trino. Consistent with previous studies of lepton number 
generation using the momentum averaged quantum rate 
equations, we find that the process is chaotic and highly 
sensitive to the mixing parameters coupling the active 
and sterile neutrino states. 

Our main results are the fractal nature of the depen- 
dence on the sign of the lepton number to the mixing 
parameters as well as the initial baryon number. This 
leaves significant implications for the predictability of the 
lepton number of the universe, even in the case of arbi- 
trarily precise determinations of the mixing parameters 



between the active and sterile states in certain regions of 
the parameter space. Future searches for subdominant 
mixings with the active neutrinos and unitarity in the 
active neutrino mixing matrix therefore could be inde- 
terminant in predicting the sign and magnitude of the 
cosmological lepton number. 

The fractal or correlation dimension depends on the 
region in parameter space that is examined, which is ex- 
pected. The banded structure of the sign of the param- 
eter space undergoes changes in its structure through- 
out the space, with some definitive lines in the pattern 
change. The fractal dimension therefore shifts in the pa- 
rameter space. This is likely due to the fact that the fre- 
quency of the rapid oscillations are deterimined by the 
oscillation dynamics, and the point at which they stop 
or freeze-out is determined by the expansion dynamics, 
which depend on the total radiation energy density of the 
universe independent of the neutrino physics. 

Our work has focused on the single-momentum quan- 
tum rate equations for the solution of this system, and 
we acknowledge that it has been shown that the behavior 
changes in the full momentum-dependent quantum ki- 
netic description. Studies have found that the oscillatory 
behavior persists and may be chaotic in the full quantum 
kinetic equation solution |23l 125] . Since there are numeri- 
cal difficulties in quantifying the nature of the solution in 
the full momentum case, which we find in our own work, 
the resolution of this question may require a completely 
new formulation of the quantum kinetic equations and 
their solution in order to make them numerically stable 
to provide a convincing solution of the evolution. The 
work presented here highlights the new type of chaotic, 
fractal physics that can arise in this system, and its con- 
sequences in determining the lepton number of the uni- 
verse given the existence of a light sterile neutrino. If a 
light sterile neutrino is discovered within this parameter 
space, a new approaches to quantifying the lepton num- 
ber of the universe will be strongly motivated. Though 
if the sign of the lepton produces such fractal structure, 
the sign of the lepton asymmetry may be indeterminate 
due to arbitrarily fine fractal structures. 
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